function [residual, g1, g2, g3] = dynaremodel_KurmannSims_Part3_dynamic(y, x, params, steady_state, it_)
%
% Status : Computes dynamic model for Dynare
%
% Inputs :
%   y         [#dynamic variables by 1] double    vector of endogenous variables in the order stored
%                                                 in M_.lead_lag_incidence; see the Manual
%   x         [nperiods by M_.exo_nbr] double     matrix of exogenous variables (in declaration order)
%                                                 for all simulation periods
%   steady_state  [M_.endo_nbr by 1] double       vector of steady state values
%   params    [M_.param_nbr by 1] double          vector of parameter values in declaration order
%   it_       scalar double                       time period for exogenous variables for which to evaluate the model
%
% Outputs:
%   residual  [M_.endo_nbr by 1] double    vector of residuals of the dynamic model equations in order of 
%                                          declaration of the equations.
%                                          Dynare may prepend auxiliary equations, see M_.aux_vars
%   g1        [M_.endo_nbr by #dynamic variables] double    Jacobian matrix of the dynamic model equations;
%                                                           rows: equations in order of declaration
%                                                           columns: variables in order stored in M_.lead_lag_incidence followed by the ones in M_.exo_names
%   g2        [M_.endo_nbr by (#dynamic variables)^2] double   Hessian matrix of the dynamic model equations;
%                                                              rows: equations in order of declaration
%                                                              columns: variables in order stored in M_.lead_lag_incidence followed by the ones in M_.exo_names
%   g3        [M_.endo_nbr by (#dynamic variables)^3] double   Third order derivative matrix of the dynamic model equations;
%                                                              rows: equations in order of declaration
%                                                              columns: variables in order stored in M_.lead_lag_incidence followed by the ones in M_.exo_names
%
%
% Warning : this file is generated automatically by Dynare
%           from model file (.mod)

%
% Model equations
%

residual = zeros(52, 1);
T19 = exp(y(70))*params(15)/params(16)*exp(y(35))^params(16);
T34 = exp(y(70))*params(15)*exp(y(35))^(params(16)-1);
T49 = y(58)^(-1);
T67 = params(12)*params(3)*exp(y(76))/(y(89)*exp(y(77))-exp(y(28))*params(12));
T73 = y(89)^(-1);
T80 = exp(y(78))^(-1);
T81 = params(3)*(1+y(29))*T73*exp(y(74))*T80;
T115 = exp(y(38))/exp(y(11));
T116 = y(58)*T115;
T120 = exp(y(25))*params(21)*(T116-params(26));
T132 = exp(y(45))/exp(y(13));
T133 = T132-1;
T136 = exp(y(36))*params(22)/2*T133^2;
T143 = exp(y(25))*(exp(y(35))*exp(y(34))*exp(y(36))-T136-exp(y(45))*exp(y(36))*params(22)*T133/exp(y(13)));
T151 = exp(y(86))/exp(y(45));
T153 = params(22)*params(3)*exp(y(74))*exp(y(81))*(T151-1);
T154 = T151^2;
T155 = T153*T154;
T167 = params(20)/2;
T171 = params(21)/2;
T176 = exp(y(82))/exp(y(39));
T177 = y(89)*T176;
T184 = exp(y(74))*(exp(y(79))*exp(y(80))-params(19)*(exp(y(80))-1)-T167*(exp(y(80))-1)^2-T171*(T177-params(26))+(T177-params(26))*params(21)*T177);
T205 = exp(y(43))^params(10);
T215 = params(5)*params(11)^((-params(10))-1);
T218 = T215*y(89)^(1+params(10));
T220 = T218*y(83);
T222 = exp(y(30));
T225 = T222^((-params(10))*params(7));
T226 = T220*T225;
T227 = exp(y(78))^params(10);
T237 = params(11)^(-params(10));
T239 = y(89)^params(10);
T244 = T222^(params(7)*(1-params(10)));
T245 = y(83)*params(5)*T237*T239*T244;
T247 = exp(y(78))^(params(10)-1);
T257 = exp(y(41))*params(10)/(params(10)-1)/exp(y(42));
T272 = T222^((-params(9))*params(8));
T273 = y(83)*y(89)*params(6)*T272;
T274 = exp(y(78))^params(9);
T283 = T222^(params(8)*(1-params(9)));
T284 = y(83)*y(89)*params(6)*T283;
T286 = exp(y(78))^(params(9)-1);
T302 = exp(y(53))^params(2);
T305 = exp(y(44))^(1-params(2));
T314 = exp(y(52))*exp(y(49))*params(2)*exp(y(53))^(params(2)-1);
T321 = exp(y(44))^(-params(2));
T322 = T302*exp(y(52))*exp(y(49))*(1-params(2))*T321;
T333 = params(6)*exp(y(5))^(params(9)*(-params(8)));
T334 = T222^params(9);
T343 = params(6)*exp(y(5))^(params(8)*(1-params(9)));
T344 = T222^(params(9)-1);
T361 = exp(y(5))^(params(10)*(-params(7)));
T362 = T237*params(5)*T239*T361;
T363 = T222^params(10);
T367 = exp(y(43))/exp(y(12));
T368 = T367^params(10);
T379 = params(11)^(1-params(10));
T380 = params(5)*y(89)^(params(10)-1)*T379;
T381 = exp(y(5))^(params(7)*(1-params(10)));
T382 = T380*T381;
T383 = T222^(params(10)-1);
T385 = exp(y(12))^(1-params(10));
T400 = exp(y(11))*T171*(T116-params(26))^2;
T407 = exp(y(11))*(params(19)*(exp(y(32))-1)+T167*(exp(y(32))-1)^2);
T409 = exp(y(58))^(-1);
T410 = T407*T409;
lhs =y(22);
rhs =params(14)+T19+params(17)/params(18)*exp(y(34))^params(18);
residual(1)= lhs-rhs;
lhs =y(23);
rhs =T34;
residual(2)= lhs-rhs;
lhs =y(24);
rhs =params(17)*exp(y(34))^(params(18)-1);
residual(3)= lhs-rhs;
lhs =exp(y(25));
rhs =exp(y(27))/(exp(y(28))-params(12)*T49*exp(y(3)))-T67;
residual(4)= lhs-rhs;
lhs =exp(y(25));
rhs =T81;
residual(5)= lhs-rhs;
lhs =exp(y(31));
rhs =params(19)+params(20)*(exp(y(32))-1);
residual(6)= lhs-rhs;
lhs =y(24)*exp(y(27))*exp(y(33));
rhs =exp(y(35))*exp(y(25))*exp(y(36));
residual(7)= lhs-rhs;
lhs =y(23)*exp(y(27))*exp(y(33));
rhs =exp(y(34))*exp(y(25))*exp(y(36));
residual(8)= lhs-rhs;
lhs =exp(y(25));
rhs =exp(y(26))*exp(y(37))-T120;
residual(9)= lhs-rhs;
lhs =y(22)*exp(y(27))*exp(y(33));
rhs =T143+T155;
residual(10)= lhs-rhs;
lhs =exp(y(26));
rhs =params(3)*T73*(T184+exp(y(75))*(1-params(4)));
residual(11)= lhs-rhs;
lhs =y(40);
rhs =exp(y(25))*T49*params(3)/exp(y(1));
residual(12)= lhs-rhs;
lhs =exp(y(41));
rhs =exp(y(36))*T205*exp(y(44))+T226*T227*exp(y(84));
residual(13)= lhs-rhs;
lhs =exp(y(42));
rhs =T205*exp(y(44))+T245*T247*exp(y(85));
residual(14)= lhs-rhs;
lhs =y(46);
rhs =T257;
residual(15)= lhs-rhs;
lhs =y(47);
rhs =exp(y(49))*exp(y(50))+T273*T274*y(87);
residual(16)= lhs-rhs;
lhs =y(48);
rhs =exp(y(50))+T284*T286*y(88);
residual(17)= lhs-rhs;
lhs =y(51);
rhs =y(47)*params(9)/(params(9)-1)/y(48);
residual(18)= lhs-rhs;
lhs =exp(y(50));
rhs =exp(y(52))*T302*T305-params(62);
residual(19)= lhs-rhs;
lhs =exp(y(31));
rhs =T305*T314;
residual(20)= lhs-rhs;
lhs =exp(y(43));
rhs =T322;
residual(21)= lhs-rhs;
lhs =y(54);
rhs =(1-params(6))*y(51)^(-params(9))+T333*T334*y(16);
residual(22)= lhs-rhs;
lhs =1;
rhs =(1-params(6))*y(51)^(1-params(9))+T343*T344;
residual(23)= lhs-rhs;
lhs =exp(y(55));
rhs =exp(y(44))*y(56);
residual(24)= lhs-rhs;
lhs =y(56);
rhs =(1-params(5))*(y(46)/exp(y(43)))^(-params(10))+T362*T363*T368*y(17);
residual(25)= lhs-rhs;
lhs =exp(y(43))^(1-params(10));
rhs =(1-params(5))*y(46)^(1-params(10))+T382*T383*T385;
residual(26)= lhs-rhs;
lhs =exp(y(55));
rhs =exp(y(45))*exp(y(35))*exp(y(34));
residual(27)= lhs-rhs;
lhs =exp(y(53));
rhs =T49*exp(y(32))*exp(y(11));
residual(28)= lhs-rhs;
lhs =exp(y(50));
rhs =exp(y(28))+exp(y(38))+exp(y(45))*T136+T49*T400+T410;
residual(29)= lhs-rhs;
lhs =exp(y(39));
rhs =exp(y(37))*exp(y(38))+T49*exp(y(11))*(1-params(4));
residual(30)= lhs-rhs;
lhs =y(29);
rhs =(1-params(28))*params(51)+params(28)*y(4)+(1-params(28))*(params(29)*(y(30)-log(params(13)))+params(30)*(y(50)-y(14)+log(y(58))-log(params(11))))+params(31)*x(it_, 2);
residual(31)= lhs-rhs;
lhs =y(57);
rhs =(1-params(33))*log(params(27))+params(33)*y(18)+params(38)*y(21);
residual(32)= lhs-rhs;
lhs =y(58);
rhs =exp(y(57))^(1/(1-params(2)));
residual(33)= lhs-rhs;
lhs =y(52);
rhs =params(32)*y(15)+params(37)*x(it_, 1);
residual(34)= lhs-rhs;
lhs =y(33);
rhs =(1-params(34))*log(params(23))+params(34)*y(7)+params(39)*x(it_, 4);
residual(35)= lhs-rhs;
lhs =y(37);
rhs =params(35)*y(10)+params(40)*x(it_, 5);
residual(36)= lhs-rhs;
lhs =y(27);
rhs =params(36)*y(2)+params(41)*x(it_, 6);
residual(37)= lhs-rhs;
lhs =y(59);
rhs =params(2)*(y(32)-y(6))+(1-params(2))*(y(34)-y(8));
residual(38)= lhs-rhs;
lhs =y(60);
rhs =y(35)-y(9);
residual(39)= lhs-rhs;
lhs =y(61);
rhs =y(60)*params(59);
residual(40)= lhs-rhs;
lhs =y(72);
rhs =exp(y(43))*exp(y(44))/exp(y(50));
residual(41)= lhs-rhs;
lhs =y(62);
rhs =y(65)-(1-y(72))*y(19)-y(72)*y(68);
residual(42)= lhs-rhs;
lhs =y(63);
rhs =y(62)-y(59);
residual(43)= lhs-rhs;
lhs =y(64);
rhs =y(62)-y(61);
residual(44)= lhs-rhs;
lhs =y(65);
rhs =y(50)-y(14)+log(y(58));
residual(45)= lhs-rhs;
lhs =y(66);
rhs =log(y(58))+y(28)-y(3);
residual(46)= lhs-rhs;
lhs =y(67);
rhs =log(y(58))+y(39)-y(11);
residual(47)= lhs-rhs;
lhs =y(68);
rhs =y(35)+y(45)-(y(13)+y(9));
residual(48)= lhs-rhs;
lhs =y(69);
rhs =y(35)+y(45);
residual(49)= lhs-rhs;
lhs =y(70);
rhs =params(61)*y(20)+params(60)*x(it_, 7);
residual(50)= lhs-rhs;
lhs =y(71);
rhs =y(52)+y(57)-y(15);
residual(51)= lhs-rhs;
lhs =y(73);
rhs =x(it_, 3);
residual(52)= lhs-rhs;
if nargout >= 2,
  g1 = zeros(52, 96);

  %
  % Jacobian matrix
  %

T596 = (exp(y(28))-params(12)*T49*exp(y(3)))*(exp(y(28))-params(12)*T49*exp(y(3)));
T605 = (y(89)*exp(y(77))-exp(y(28))*params(12))*(y(89)*exp(y(77))-exp(y(28))*params(12));
T782 = y(58)*(-(exp(y(38))*exp(y(11))))/(exp(y(11))*exp(y(11)));
T798 = y(89)*(-(exp(y(82))*exp(y(39))))/(exp(y(39))*exp(y(39)));
T837 = getPowerDeriv(T367,params(10),1);
T847 = exp(y(43))*getPowerDeriv(exp(y(43)),params(10),1);
T857 = getPowerDeriv(y(46)/exp(y(43)),(-params(10)),1);
T870 = exp(y(44))*getPowerDeriv(exp(y(44)),1-params(2),1);
T883 = (-(exp(y(45))*exp(y(13))))/(exp(y(13))*exp(y(13)));
T915 = (-(exp(y(45))*exp(y(86))))/(exp(y(45))*exp(y(45)));
T971 = exp(y(53))*getPowerDeriv(exp(y(53)),params(2),1);
T990 = getPowerDeriv(y(58),(-1),1);
T1029 = getPowerDeriv(y(89),(-1),1);
T1054 = getPowerDeriv(y(89),params(10),1);
  g1(1,22)=1;
  g1(1,34)=(-(params(17)/params(18)*exp(y(34))*getPowerDeriv(exp(y(34)),params(18),1)));
  g1(1,35)=(-(exp(y(70))*params(15)/params(16)*exp(y(35))*getPowerDeriv(exp(y(35)),params(16),1)));
  g1(1,70)=(-T19);
  g1(2,23)=1;
  g1(2,35)=(-(exp(y(70))*params(15)*exp(y(35))*getPowerDeriv(exp(y(35)),params(16)-1,1)));
  g1(2,70)=(-T34);
  g1(3,24)=1;
  g1(3,34)=(-(params(17)*exp(y(34))*getPowerDeriv(exp(y(34)),params(18)-1,1)));
  g1(4,25)=exp(y(25));
  g1(4,27)=(-(exp(y(27))/(exp(y(28))-params(12)*T49*exp(y(3)))));
  g1(4,76)=T67;
  g1(4,3)=(-((-(exp(y(27))*(-(params(12)*T49*exp(y(3))))))/T596));
  g1(4,28)=(-((-(exp(y(27))*exp(y(28))))/T596-(-(params(12)*params(3)*exp(y(76))*(-(exp(y(28))*params(12)))))/T605));
  g1(4,77)=(-(params(12)*params(3)*exp(y(76))*y(89)*exp(y(77))))/T605;
  g1(4,58)=(-((-(exp(y(27))*(-(exp(y(3))*params(12)*T990))))/T596));
  g1(4,89)=(-(params(12)*params(3)*exp(y(76))*exp(y(77))))/T605;
  g1(5,25)=exp(y(25));
  g1(5,74)=(-T81);
  g1(5,29)=(-(T80*exp(y(74))*params(3)*T73));
  g1(5,78)=(-(params(3)*(1+y(29))*T73*exp(y(74))*exp(y(78))*getPowerDeriv(exp(y(78)),(-1),1)));
  g1(5,89)=(-(T80*exp(y(74))*params(3)*(1+y(29))*T1029));
  g1(6,31)=exp(y(31));
  g1(6,32)=(-(params(20)*exp(y(32))));
  g1(7,24)=exp(y(27))*exp(y(33));
  g1(7,25)=(-(exp(y(35))*exp(y(25))*exp(y(36))));
  g1(7,27)=y(24)*exp(y(27))*exp(y(33));
  g1(7,33)=y(24)*exp(y(27))*exp(y(33));
  g1(7,35)=(-(exp(y(35))*exp(y(25))*exp(y(36))));
  g1(7,36)=(-(exp(y(35))*exp(y(25))*exp(y(36))));
  g1(8,23)=exp(y(27))*exp(y(33));
  g1(8,25)=(-(exp(y(34))*exp(y(25))*exp(y(36))));
  g1(8,27)=y(23)*exp(y(27))*exp(y(33));
  g1(8,33)=y(23)*exp(y(27))*exp(y(33));
  g1(8,34)=(-(exp(y(34))*exp(y(25))*exp(y(36))));
  g1(8,36)=(-(exp(y(34))*exp(y(25))*exp(y(36))));
  g1(9,25)=exp(y(25))-(-T120);
  g1(9,26)=(-(exp(y(26))*exp(y(37))));
  g1(9,37)=(-(exp(y(26))*exp(y(37))));
  g1(9,38)=exp(y(25))*params(21)*T116;
  g1(9,11)=exp(y(25))*params(21)*T782;
  g1(9,58)=exp(y(25))*params(21)*T115;
  g1(10,22)=exp(y(27))*exp(y(33));
  g1(10,25)=(-T143);
  g1(10,74)=(-T155);
  g1(10,27)=y(22)*exp(y(27))*exp(y(33));
  g1(10,33)=y(22)*exp(y(27))*exp(y(33));
  g1(10,34)=(-(exp(y(25))*exp(y(35))*exp(y(34))*exp(y(36))));
  g1(10,35)=(-(exp(y(25))*exp(y(35))*exp(y(34))*exp(y(36))));
  g1(10,36)=(-T143);
  g1(10,81)=(-T155);
  g1(10,13)=(-(exp(y(25))*((-(exp(y(36))*params(22)/2*T883*2*T133))-(exp(y(13))*exp(y(45))*exp(y(36))*params(22)*T883-exp(y(13))*exp(y(45))*exp(y(36))*params(22)*T133)/(exp(y(13))*exp(y(13))))));
  g1(10,45)=(-(exp(y(25))*((-(exp(y(36))*params(22)/2*T132*2*T133))-(exp(y(45))*exp(y(36))*params(22)*T133+exp(y(45))*exp(y(36))*params(22)*T132)/exp(y(13)))+T154*params(22)*params(3)*exp(y(74))*exp(y(81))*T915+T153*T915*2*T151));
  g1(10,86)=(-(T154*params(22)*params(3)*exp(y(74))*exp(y(81))*T151+T153*T151*2*T151));
  g1(11,74)=(-(params(3)*T73*T184));
  g1(11,26)=exp(y(26));
  g1(11,75)=(-(params(3)*T73*exp(y(75))*(1-params(4))));
  g1(11,79)=(-(params(3)*T73*exp(y(74))*exp(y(79))*exp(y(80))));
  g1(11,80)=(-(params(3)*T73*exp(y(74))*(exp(y(79))*exp(y(80))-params(19)*exp(y(80))-T167*exp(y(80))*2*(exp(y(80))-1))));
  g1(11,82)=(-(params(3)*T73*exp(y(74))*((-(T171*T177))+(T177-params(26))*params(21)*T177+T177*params(21)*T177)));
  g1(11,39)=(-(params(3)*T73*exp(y(74))*((-(T171*T798))+params(21)*T177*T798+(T177-params(26))*params(21)*T798)));
  g1(11,89)=(-((T184+exp(y(75))*(1-params(4)))*params(3)*T1029+params(3)*T73*exp(y(74))*((-(T171*T176))+T176*params(21)*T177+(T177-params(26))*params(21)*T176)));
  g1(12,1)=(-((-(exp(y(25))*T49*params(3)*exp(y(1))))/(exp(y(1))*exp(y(1)))));
  g1(12,25)=(-(exp(y(25))*T49*params(3)/exp(y(1))));
  g1(12,40)=1;
  g1(12,58)=(-(exp(y(25))*params(3)*T990/exp(y(1))));
  g1(13,30)=(-(exp(y(84))*T227*T220*T222*getPowerDeriv(T222,(-params(10))*params(7),1)));
  g1(13,78)=(-(exp(y(84))*T226*exp(y(78))*getPowerDeriv(exp(y(78)),params(10),1)));
  g1(13,36)=(-(exp(y(36))*T205*exp(y(44))));
  g1(13,83)=(-(exp(y(84))*T227*T218*T225));
  g1(13,41)=exp(y(41));
  g1(13,84)=(-(T226*T227*exp(y(84))));
  g1(13,43)=(-(exp(y(44))*exp(y(36))*T847));
  g1(13,44)=(-(exp(y(36))*T205*exp(y(44))));
  g1(13,89)=(-(exp(y(84))*T227*T225*y(83)*T215*getPowerDeriv(y(89),1+params(10),1)));
  g1(14,30)=(-(exp(y(85))*T247*y(83)*params(5)*T237*T239*T222*getPowerDeriv(T222,params(7)*(1-params(10)),1)));
  g1(14,78)=(-(exp(y(85))*T245*exp(y(78))*getPowerDeriv(exp(y(78)),params(10)-1,1)));
  g1(14,83)=(-(exp(y(85))*T247*params(5)*T237*T239*T244));
  g1(14,42)=exp(y(42));
  g1(14,85)=(-(T245*T247*exp(y(85))));
  g1(14,43)=(-(exp(y(44))*T847));
  g1(14,44)=(-(T205*exp(y(44))));
  g1(14,89)=(-(exp(y(85))*T247*T244*y(83)*params(5)*T237*T1054));
  g1(15,41)=(-T257);
  g1(15,42)=(-((-(exp(y(42))*exp(y(41))*params(10)/(params(10)-1)))/(exp(y(42))*exp(y(42)))));
  g1(15,46)=1;
  g1(16,30)=(-(y(87)*T274*y(83)*y(89)*params(6)*T222*getPowerDeriv(T222,(-params(9))*params(8),1)));
  g1(16,78)=(-(y(87)*T273*exp(y(78))*getPowerDeriv(exp(y(78)),params(9),1)));
  g1(16,83)=(-(y(87)*T274*y(89)*params(6)*T272));
  g1(16,47)=1;
  g1(16,87)=(-(T273*T274));
  g1(16,49)=(-(exp(y(49))*exp(y(50))));
  g1(16,50)=(-(exp(y(49))*exp(y(50))));
  g1(16,89)=(-(y(87)*T274*T272*y(83)*params(6)));
  g1(17,30)=(-(y(88)*T286*y(83)*y(89)*params(6)*T222*getPowerDeriv(T222,params(8)*(1-params(9)),1)));
  g1(17,78)=(-(y(88)*T284*exp(y(78))*getPowerDeriv(exp(y(78)),params(9)-1,1)));
  g1(17,83)=(-(y(88)*T286*y(89)*params(6)*T283));
  g1(17,48)=1;
  g1(17,88)=(-(T284*T286));
  g1(17,50)=(-exp(y(50)));
  g1(17,89)=(-(y(88)*T286*T283*y(83)*params(6)));
  g1(18,47)=(-(params(9)/(params(9)-1)/y(48)));
  g1(18,48)=(-((-(y(47)*params(9)/(params(9)-1)))/(y(48)*y(48))));
  g1(18,51)=1;
  g1(19,44)=(-(exp(y(52))*T302*T870));
  g1(19,50)=exp(y(50));
  g1(19,52)=(-(exp(y(52))*T302*T305));
  g1(19,53)=(-(T305*exp(y(52))*T971));
  g1(20,31)=exp(y(31));
  g1(20,44)=(-(T314*T870));
  g1(20,49)=(-(T305*T314));
  g1(20,52)=(-(T305*T314));
  g1(20,53)=(-(T305*exp(y(52))*exp(y(49))*params(2)*exp(y(53))*getPowerDeriv(exp(y(53)),params(2)-1,1)));
  g1(21,43)=exp(y(43));
  g1(21,44)=(-(T302*exp(y(52))*exp(y(49))*(1-params(2))*exp(y(44))*getPowerDeriv(exp(y(44)),(-params(2)),1)));
  g1(21,49)=(-T322);
  g1(21,52)=(-T322);
  g1(21,53)=(-(T321*exp(y(52))*exp(y(49))*(1-params(2))*T971));
  g1(22,5)=(-(y(16)*T334*params(6)*exp(y(5))*getPowerDeriv(exp(y(5)),params(9)*(-params(8)),1)));
  g1(22,30)=(-(y(16)*T333*T222*getPowerDeriv(T222,params(9),1)));
  g1(22,51)=(-((1-params(6))*getPowerDeriv(y(51),(-params(9)),1)));
  g1(22,16)=(-(T333*T334));
  g1(22,54)=1;
  g1(23,5)=(-(T344*params(6)*exp(y(5))*getPowerDeriv(exp(y(5)),params(8)*(1-params(9)),1)));
  g1(23,30)=(-(T343*T222*getPowerDeriv(T222,params(9)-1,1)));
  g1(23,51)=(-((1-params(6))*getPowerDeriv(y(51),1-params(9),1)));
  g1(24,44)=(-(exp(y(44))*y(56)));
  g1(24,55)=exp(y(55));
  g1(24,56)=(-exp(y(44)));
  g1(25,5)=(-(y(17)*T368*T363*T237*params(5)*T239*exp(y(5))*getPowerDeriv(exp(y(5)),params(10)*(-params(7)),1)));
  g1(25,30)=(-(y(17)*T368*T362*T222*getPowerDeriv(T222,params(10),1)));
  g1(25,12)=(-(y(17)*T362*T363*(-(exp(y(43))*exp(y(12))))/(exp(y(12))*exp(y(12)))*T837));
  g1(25,43)=(-((1-params(5))*(-(exp(y(43))*y(46)))/(exp(y(43))*exp(y(43)))*T857+y(17)*T362*T363*T367*T837));
  g1(25,46)=(-((1-params(5))*T857*1/exp(y(43))));
  g1(25,17)=(-(T362*T363*T368));
  g1(25,56)=1;
  g1(25,89)=(-(y(17)*T368*T363*T361*T237*params(5)*T1054));
  g1(26,5)=(-(T385*T383*T380*exp(y(5))*getPowerDeriv(exp(y(5)),params(7)*(1-params(10)),1)));
  g1(26,30)=(-(T385*T382*T222*getPowerDeriv(T222,params(10)-1,1)));
  g1(26,12)=(-(T382*T383*exp(y(12))*getPowerDeriv(exp(y(12)),1-params(10),1)));
  g1(26,43)=exp(y(43))*getPowerDeriv(exp(y(43)),1-params(10),1);
  g1(26,46)=(-((1-params(5))*getPowerDeriv(y(46),1-params(10),1)));
  g1(26,89)=(-(T385*T383*T381*T379*params(5)*getPowerDeriv(y(89),params(10)-1,1)));
  g1(27,34)=(-(exp(y(45))*exp(y(35))*exp(y(34))));
  g1(27,35)=(-(exp(y(45))*exp(y(35))*exp(y(34))));
  g1(27,45)=(-(exp(y(45))*exp(y(35))*exp(y(34))));
  g1(27,55)=exp(y(55));
  g1(28,32)=(-(T49*exp(y(32))*exp(y(11))));
  g1(28,11)=(-(T49*exp(y(32))*exp(y(11))));
  g1(28,53)=exp(y(53));
  g1(28,58)=(-(exp(y(32))*exp(y(11))*T990));
  g1(29,28)=(-exp(y(28)));
  g1(29,32)=(-(T409*exp(y(11))*(params(19)*exp(y(32))+T167*exp(y(32))*2*(exp(y(32))-1))));
  g1(29,36)=(-(exp(y(45))*T136));
  g1(29,38)=(-(exp(y(38))+T49*exp(y(11))*T171*T116*2*(T116-params(26))));
  g1(29,11)=(-(T410+T49*(T400+exp(y(11))*T171*2*(T116-params(26))*T782)));
  g1(29,13)=(-(exp(y(45))*exp(y(36))*params(22)/2*T883*2*T133));
  g1(29,45)=(-(exp(y(45))*T136+exp(y(45))*exp(y(36))*params(22)/2*T132*2*T133));
  g1(29,50)=exp(y(50));
  g1(29,58)=(-(T400*T990+T49*exp(y(11))*T171*T115*2*(T116-params(26))+T407*exp(y(58))*getPowerDeriv(exp(y(58)),(-1),1)));
  g1(30,37)=(-(exp(y(37))*exp(y(38))));
  g1(30,38)=(-(exp(y(37))*exp(y(38))));
  g1(30,11)=(-(T49*exp(y(11))*(1-params(4))));
  g1(30,39)=exp(y(39));
  g1(30,58)=(-(exp(y(11))*(1-params(4))*T990));
  g1(31,4)=(-params(28));
  g1(31,29)=1;
  g1(31,30)=(-((1-params(28))*params(29)));
  g1(31,14)=(-((1-params(28))*(-params(30))));
  g1(31,50)=(-((1-params(28))*params(30)));
  g1(31,58)=(-((1-params(28))*params(30)*1/y(58)));
  g1(31,91)=(-params(31));
  g1(32,18)=(-params(33));
  g1(32,57)=1;
  g1(32,21)=(-params(38));
  g1(33,57)=(-(exp(y(57))*getPowerDeriv(exp(y(57)),1/(1-params(2)),1)));
  g1(33,58)=1;
  g1(34,15)=(-params(32));
  g1(34,52)=1;
  g1(34,90)=(-params(37));
  g1(35,7)=(-params(34));
  g1(35,33)=1;
  g1(35,93)=(-params(39));
  g1(36,10)=(-params(35));
  g1(36,37)=1;
  g1(36,94)=(-params(40));
  g1(37,2)=(-params(36));
  g1(37,27)=1;
  g1(37,95)=(-params(41));
  g1(38,6)=params(2);
  g1(38,32)=(-params(2));
  g1(38,8)=1-params(2);
  g1(38,34)=(-(1-params(2)));
  g1(38,59)=1;
  g1(39,9)=1;
  g1(39,35)=(-1);
  g1(39,60)=1;
  g1(40,60)=(-params(59));
  g1(40,61)=1;
  g1(41,43)=(-(exp(y(43))*exp(y(44))/exp(y(50))));
  g1(41,44)=(-(exp(y(43))*exp(y(44))/exp(y(50))));
  g1(41,50)=(-((-(exp(y(50))*exp(y(43))*exp(y(44))))/(exp(y(50))*exp(y(50)))));
  g1(41,72)=1;
  g1(42,62)=1;
  g1(42,65)=(-1);
  g1(42,19)=1-y(72);
  g1(42,68)=y(72);
  g1(42,72)=(-(y(19)-y(68)));
  g1(43,59)=1;
  g1(43,62)=(-1);
  g1(43,63)=1;
  g1(44,61)=1;
  g1(44,62)=(-1);
  g1(44,64)=1;
  g1(45,14)=1;
  g1(45,50)=(-1);
  g1(45,58)=(-(1/y(58)));
  g1(45,65)=1;
  g1(46,3)=1;
  g1(46,28)=(-1);
  g1(46,58)=(-(1/y(58)));
  g1(46,66)=1;
  g1(47,11)=1;
  g1(47,39)=(-1);
  g1(47,58)=(-(1/y(58)));
  g1(47,67)=1;
  g1(48,9)=1;
  g1(48,35)=(-1);
  g1(48,13)=1;
  g1(48,45)=(-1);
  g1(48,68)=1;
  g1(49,35)=(-1);
  g1(49,45)=(-1);
  g1(49,69)=1;
  g1(50,20)=(-params(61));
  g1(50,70)=1;
  g1(50,96)=(-params(60));
  g1(51,15)=1;
  g1(51,52)=(-1);
  g1(51,57)=(-1);
  g1(51,71)=1;
  g1(52,92)=(-1);
  g1(52,73)=1;

if nargout >= 3,
  %
  % Hessian matrix
  %

  g2 = sparse([],[],[],52,9216);
if nargout >= 4,
  %
  % Third order derivatives
  %

  g3 = sparse([],[],[],52,884736);
end
end
end
end
